#if !defined(WRAPPER_QUASI)
#define WRAPPER_QUASI

#include <iostream>
#include "function.h"
#include "nr_cpp/nr3.h"
#include "nr_cpp/ludcmp.h"
#include "nr_cpp/qrdcmp.h"
#include "nr_cpp/roots_multidim.h"
#include "nr_cpp/quasinewton.h"
#include <string>

struct Funcdd : public Funcd{

  Funcdd(int numdim, Function *func)
  {
    NInput = numdim;
    function = func;

    newA(x_new, NInput);
    newA(deriv_new, NInput);
  }
  ~Funcdd(){}
  
  // print out solution to stream
  void print_sol(VecDoub_I &x, double ChiSq, ofstream& append, string cifLabel = "default")
  {
    for (int i = 0; i < NInput; i++)
    {
      x_new[i+1] = x[i];
    }
    function->print_sol(x_new, ChiSq, append, cifLabel);
  }

  Doub operator() ( VecDoub_I &x )
  {
    for (int i = 0; i < NInput; i++)
    {
      x_new[i+1] = x[i];
    }

    return function->Value(x_new);
  }

  void df(VecDoub_I &x, VecDoub_O &deriv)
  {
    for (int i = 0; i < NInput; i++)
    {
      x_new[i+1] = x[i];
      deriv_new[i+1] = deriv[i];
    }

    function->BuildDerivative(x_new, deriv_new);

    for (int i = 0; i < NInput; i++)
    {
      deriv[i] = deriv_new[i+1];
    }
  }

  int NInput;

private:
  double *x_new;
  double *deriv_new;
  Function *function;
};

/** Wrapper for Numerical recipies dfpmin function (QN) 
 *
 * @param pI InOut parameter. As input: starting point for minimization and as output: final point of mimimization 
 * @param gtol Input parameter specifies convergence criterion on when gradient close enough to zero
 * @param fret Output value of final minimization function value
 * @param Funcdd Function to minimize
 *  
 */
void WrapperQuasi(double *pI, double gtol, double &fret, Funcdd &funcd, string& directory);

#endif // !defined(WRAPPER_QUASI)